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ABSTRACT 

We present results from the first numerical analysis to support the hypothe- 
sis, first proposed in Coleman et al. (2004), that the Fornax dwarf galaxy was 
formed from the minor merging of two dwarfs about 2 Gyr ago. Using orbits 
\ for the Fornax dwarf that are consistent with the latest proper motion mea- 

^ surements, our dynamical evolution models show that the observed asymmetric 

shell-like substructures can be formed from the remnant of a smaller dwarf dur- 
ing minor merging. These models also predict the formation of diffuse stellar 
streams. We discuss how these stellar substructures depend on model parame- 
ters of dwarf-dwarf merging, and how the intermediate-age subpopulations found 
in the vicinity of these substructures may be formed from gas accretion in the 
past merger events. We also suggest that one of Fornax's GCs originates from a 
merged dwarf companion, and demonstrate where as yet undetected tidal streams 
or HI gas formed from the dwarf merging may be found in the outer halo of the 
Galaxy. 

Subject headings: galaxies: kinematics and dynamics — galaxies: individual 
(Fornax) — Local Group 



Introduction 



The Fornax galaxy possesses many unusual characteristics that have warranted far 
greater investigation than might be expected for a typical dwarf spheroidal (dSph). Lying 
in the Local Group at a distance of 138kpc (Saviane et al. 2000) and receding at 53.3kms~^ 
(Walker et al. 2007), the foremost distinguishing features are six Globular Clusters (GC), 
resulting in one of the highest known specific frequencies of all galaxies (29, from van den 
Bergh 1998). Although these GCs are associated with Fornax (Mateo et al. 1991), with 
distances to the galaxy center ranging between 0.24 and 1.60kpc, there exists an atypical 
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spread in metallicity and age (Buonanno et al. 1999, and Strader et al. 2003). Moreover, 
dynamical friction acting on the GCs in the dark matter (DM) dominated Fornax should 
have resulted in their assimilation into the Fornax nucleus since their formation more than 
~10 Gyr ago. Tidal interaction with the Milky Way (MW) has been invoked to explain the 
latter by Oh et al. (2000); further studies by Goerdt et al. (2006) propose that a cored DM 
halo would increase the effective sinking time for GCs, while Cowsik et al. (2009) asserts 
dynamical friction to be weaker than previously thought. 

The above discussion lies in parallel but separate from recent photometric analyses by 
Coleman ct al. (2004, 2005). These studies have interpreted statistically-significant {4.5a) 
stellar overdensities as evidence of an extant shell substructure, with the most prominent 
shell 'sheets' aligned with the major axis of Fornax and extending beyond the tidal radius 
{vf = 71 ± 4 arcmin, Irwin & Hatzidimitriou 1995) of Fornax . Tidal interaction with the 
MW is discarded as a causative factor, since the predominantly intermediate (2 Gyr) age of 
the overdensity subpopulation contrasts with the otherwise mixed aged general population 
of Fornax (Buonanno et al. 1999, Stetson et al. 1998). Quinn (1984) showed that low-energy 
minor mergers can result in shell-like substructures, and a gas-rich companion galaxy may 
also provide the fuel source for a past starburst, given the present lack of Hj gas observed 
in Fornax. Olszewski et al. (2006) argue however that the metallicity in the overdensities 
bear the signature of Fornax itself. 

In both debates, the infiuence of the Milky Way (MW) potential has been contemplated 
with respect to the assumed orbit of Fornax. Successive measurements for the present day 
proper motion of Fornax have been conducted recently by Dinescu et al. (2004), Piatek et al. 
(2007), Walker et al. (2008), and Mendez et al. (2011), hereafter D04, P07, W08 and Mil 
respectively. Mil gives a high-energy orbit with an orbital period ~21 Gyr, while D04, P07 
and W08 imply a low-energy orbit with orbital periods ranging between 5-8 Gyr. In all cases, 
Fornax is conceived to have been at perigalacticon in the last ~ 100-300 Gyr. Mendez et al. 
(2011) propose that this triggered an enhanced star formation at this period, and perhaps a 
similarly induced starburst >10 Gyr ago, thus supporting their hyper-extended orbit (if the 
merger scenario can explain the remaining intermediate subpopulation). By contrast, the 
minimum galactocentric radial distance of 131-148kpc leads Coleman et al. (2005) to assert 
the relative insignificance of tidal forces at this separation. To date, no numerical studies 
have been conducted to assert the interaction and merger history of Fornax, using the latest 
measurements of its proper motion. 

The purpose of this LETTER is to test the hypothesis of Coleman et al. (2004) that 
Fornax is a remnant of recent dwarf merging. We show that the overdensities are evidence 
of a shell substructure formed by the most recent merger event. We also discuss how the 
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recent dwarf merging is important for (i) the formation of intermediate- age (~2 Gyr) stars 
in the core of Fornax and (ii) physical properties of its unusual globular cluster system. We 
make predictions for future observations that would constrain the present range of proper 
motions and further support the merger hypothesis. 



2. The model 

Our investigation assumes that the present day Fornax dSph is the product of the merger 
of two dwarfs, referred to hereafter as Fornax A (host) and B (companion) respectively. Our 
merger model is developed in two stages: first, we evaluate the orbital history of Fornax for 
a given Galactic potential and a mass model for Fornax A/B (§2.1) by using the "backward 
integration scheme" used in similar simulations (e.g., Diaz & Bekki 2012, hereafter DB12). 
In the second step (§2.2) we introduce an coUisionless N-body representation of Fornax B, 
and investigate its dynamical evolution under the gravitational influence of the Galaxy and 
Fornax A. 



2.1. Orbit Integration 

We adopt the same Galaxy model as used in our previous simulations on the cvohition 
of the Magellanic Cloud (DB12). The dominant component is the extended halo, for which 
we adopt a NFW density distribution (Navarro et al. 1996), up to a virial radius Ryir of 
175kpc: 

PNFwyr) 



from which the virial mass M^r is derived as: 

M^ir = AT:pQ{Ryir/c^){ln{l + c) - c/(l + c)) = 1.30 x IO^^Mq, 

where r is the spherical radius, po is the characteristic density and c is a concentration 
parameter (= 12). The disk uses the potential of Miyamoto & Nagai (1975), which has the 
form: 

Vi?2 + (a + VzTThy 
where the total mass of the disk = 5.0 x IO^^Mq (Binney & Tremaine 2008), R = 
^/ x'^ + y'^, and a and b are scale parameters with values 3.5kpc and 0.35kpc respectively. 
Finally, the Galactic bulge uses the spherical potential model of Hernquist (1990): 

GM, 

^bulge — ; 

r + Cb 
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where the mass of the bulge is 0.5 x lO^^M© and Cb — 0.7kpc (Binney & Tremaine 
2008). The total mass of the MW, up to a limit radius of 300kpc, is thus calculated as 
~ 1.73 X IO^^Mq. Fornax A is assumed to have a Plummer potential: 

d> ( ^ 



2' 



Vr2 + apA 

where MpA is the total mass of (lO^M©, consistent with the latest observations (e.g., de Boer 
et al. 2011), and apA is a softening radius of O.Gkpc. We consider the spherical Plummer 
model sufficient for the shghtly oblate extant Fornax dSph; Hernquist & Quinn (1989) show 
a complex dependence of the shell structure alignment on primary galaxy potential, but we 
consider this dependence sufficiently small to ignore within the scope of this preliminary 
study. 

The time integration of the equation of motion is performed by using 2nd-order leap-frog 
method (Murai & Fujimoto 1980) with a time step interval of ~ 1.4 x 10^ yr. For initial 
space velocities, we Tised the proper motions of P07 and Mil as representative of the low- 
and high-energy orbits respectively. The calculated error in the measurements of P07 were 
sufficiently low to justify the use of space velocities (in the galactocentric rest frame) derived 
from the nominal proper motion. The measurements of Mil by contrast are accompanied 
by significant confidence limits, in which case we separately considered the nominal (Mil), 
lowest- (MllL) and highest- (MllH) energy combinations from this data to attain the corre- 
sponding space velocities {Vx, Vy, Vz) in Table 1. We use the coordinate system adopted in 
DB12, in which Fornax is currently located at position {X,Y,Z) = (-39.0, -48.0, -126.0) 
kpc, as derived from van den Bergh (1999). For each orbit, the integration scheme is utilised 
until the time {Tm) at which a merger is judged to have commenced. To establish the time 
varying dependence of our postulated mergers, and accommodate the age range of yoimg 
stars in Fornax (e.g. Fig.4, Tsujimoto 2011), we vary the merger time among models from 
Tm = -3.5 to -2 Gyr, where T = is the present time. 



2.2. N-Body model and merger parameters 

Fornax B is assumed to be a dwarf disk galaxy embedded in a massive dark matter halo, 
and is therefore represented here by N-body particles in a self-consistent manner. We adopt 
the dwarf disk model adopted in our previous simulations on dynamical evolution of dwarf 
disk galaxies around the Galaxy, the details of which are given in Bekki & Yong (2012). 

The total mass of Fornax B (Mt^dw) is assumed as 5% of Fornax A (i.e., Mt^dw = 
5 X lO^M©). The observed average metallicity of Fornax is consistent with the merger of 
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two satellites (Fornax A and B) in a 20:1 stellar mass ratio assuming that the two merging 
satellites conform to the luminosity- metalli city relationship of MW dSph satellites (e.g., 
Kirby et al. 2011). The density profile of the dark matter halo with a total mass of Mdm.dw 
is represented by that proposed by Salucci & Burkert (2000): 

/ _ Pdm,0 

where Pdm,o and Cdm are the central dark matter density and the core (scale) radius, respec- 
tively. The mass-ratio of Mdm.dw to Mt^dw is set to be 0.9, while adm is varied for different 
models. The stellar component of the dwarf is modelled as a bulge-less stellar disk with 
the total mass of Ms^dw and the size of Rs,dw The radial (R) and vertical (Z) density pro- 
files of the stellar disk arc assumed to be proportional to exp{—R/Ro) with scale length 
Rq = 0.2i?s^dw and to sech^(Z/Zo) with scale length Zq = 0.04i?s,dwi respectively. The 
gravitational softening length for the dark matter and stellar disk is set to be 2.7pc. 

In addition to the rotational velocity caused by the gravitational field of disk and dark 
halo components, the initial radial and azimuthal velocity dispersions are assigned to the 
disc component according to the cpicyclic theory with Toomre's parameter Q = 1.5. The 
mass-ratio of Ms^dw to Mt,dw is set to be 0.1. The spin of the dwarf galaxy is specified by 
two angles 9 and (p (in units of degrees). Here, 9 is the angle between the z-axis and the 
vector of the angular momentum of the disk, and (p is the azimuthal angle measured from 
the X-axis to the projection of the angular momentum vector of the disk onto the x-y plane. 
The initial position of Fornax B, with respect to Fornax A, is described with: 

^AB = (rAscos($),0,r^i?sin($)), 

where tab is the initial radial separation and $ is an inclination angle; similarly, the velocity 
is: 

yAB^{frVcJtVc,0), 

where is the circular velocity at tab, and fr and ft are radial and tangential component 
factors respectively. The locations of the overdensities from Coleman et al. (2004, 2005) 
would appear to describe an aligned shell, which is most often imposed by radial or near- 
radial collisions (Hernquist & Quinn 1989). The relative sparscncss of the Local Group, and 
Fornax's location on the periphery of those satellites considered bound to the Galaxy, would 
result in galaxy interactions tending towards those observed in field rather than cluster 
environments. We thus constrained the trajectory of Fornax B with respect to A to a 
low-energy, infalling parabohc orbit, with fr — —0.5 and ft — 0.1. The morphology of a 
shell substructure is also sensitive to the initial orientation of the encounter: the sharp-edged 
features of the shell are phase- wrapped sheets as seen edge-on, with their overall distribution 
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(the symmetry plane) seen face-on. We assumed that the overdensities fit this perspective, 
and thus we determined ^ = 135, = and $ = 30°, which together promote a face-on view 
from a hehocentric location at T = 0. 

The dynamical evolution of Fornax B was investigated for a wide extent of the contingent 
parameter space described by the variables here defined, with models consisting of 200,000 
particles (split equally between the stellar and DM components). For parameter sets where 
the stellar distribution is comparable with extant observations, we confirm the results with 
higher resolution 1,000,000 particle models. All simulations are implemented with an original 
code (e.g., Bckki 2011), which executes with the parallel computing architecture CUDA, on 
a PC cluster installed with Nvidia GTX580 CPUs at the University of Western Australia. 

3. Results 

Figure 1 shows the backwards integrated orbits of Fornax A for the nominal proper 
motions of P07 and Mil, and conveys the wide discrepancy in their respective proximity 
to the MW. For a merger instantiated at {Tm = —3.5) using the proper motion of P07 we 
find, in the dynamic stellar distribution, overdensities are developed for a specific parameter 
subset of Fornax B that we define as Model 1 (Table 1). The top row of Figure 2 gives the 
projected stellar distribution of Model 1, with a Gaussian-convolved image of Fornax from 
Coleman et al. (2005) for comparison; the outer overdensity is highhghted together with the 
inferred axis of symmetry. At T = —1.4, a classical type 1 (Prieur 1990 classification) 'bow- 
tie' structure is prominent. The regions of maximal stellar density are further elucidated 
by first dividing the on-sky distribution projection into quadrants, and then counting the 
number of stars that lie within each equal width bin of radius from Fornax A. Accordingly, 
the top-left panel of Figure 3 shows salient maximal densities at 1° to 1.5° arranged in 
opposing quadrants, analogous to Coleman et al. (2005). As in other studies (e.g. Weil et 
al 1996), these peaks exist briefly (~0.1 Gyr) before diminishing through mixing. Fornax 
B is not completely disrupted, but we find it to first reside at the maximal stellar density 
region, before later descending to Fornax A's stellar core. Notably, we find a bimodality in 
the heliocentric velocity of the Fornax B stellar remnant up to a radius of ~ 0.7°; this can 
be compared with a similar bimodality observed for the metal-poor stars of the Fornax dSph 
(Battaglia et al. 2006), as conveyed in the right-side panels of Fig. 3). 

When an otherwise identical model is run but with Tm = 2.1 Gyr (Model 2), the 
anticipated overdensities are in fact not reproduced at T = 0. Tidal influence on Fornax is 
deemed strongest therefore during its approach to perigalacticon, and to recover overdensities 
at T = for this Tm, Fornax B is required to be more compact in its original state. With 
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a decreased disk -Rg.dw and halo core radius adm, Model 3 is conveyed in Fig. 2 and Fig. 3 
to possess a distinct overdense remnant (the surviving stellar core of Fornax B) at T = 0, 
which is again comparable qualitatively with Fornax. Taking a broader view of the stellar 
distribution, the left-side panels of Figure 4 show a tidal stream in Models 1 and 3, which 
extends up to a length of ~ 30° in alignment with the prior orbit of Fornax, and for a further 
~ 20° in the opposing direction. 

In contrast, our models with high energy orbits (Mil, MllL and MllH) develop largely 
axisymmetric shell structures (Type 2), as conveyed by the representative Model 4 (with 
proper motion Mil) of Fig. 2 (lower-right). This variation from Models 1-3 is attributed to 
the lessened capacity of the Galaxy to lock Fornax B into a more radial orbit and accretion 
by its host. These models have no tidal streams; in the case where tidal influence is most 
endorsed, with MllL and a lesser initial density of Fornax B (Model 5), small tidal plumes 
do develop (Fig. 4, lower-right panel). 

4. Discussion and Conclusions 

In this letter, we have demonstrated that the Fornax galaxy could be the remnant 
of a recent merger. We beheve that Models 1 and 3 define two bounding cases, where a 
shell substructure, qualitatively similar to that observed in the Fornax dSph, is achieved 
for a merger occurring between Tm — -3.5 and Tm — -2.1 Gyr ago. Future photometry of 
these proposed shell regions may reveal their full distribution and thus allow further fine 
tuning of the merger time (and the initial density Fornax B); at this stage, however, the 
merger hypothesis put forward by Coleman et al. (2004) appears valid. Furthermore, the 
morphology of this structure may be indicative of the orbital motion of Fornax; we anticipate 
the discovery of more substructure, but the current observations are most commensurate 
with the low-energy orbit of P07 (implying also D04, W08). We have additionally predicted 
tidal structure, similar to the well-studied stellar streams of the Sagittarius dSph (Martinez- 
Delgado et al. 2004). Future observations of extant tails can provide constraints therefore 
on the proper motion of Fornax. We suggest further that the maximal density regions of 
the substructure could represent the core of the companion, which has been not completely 
disrupted by the merger and subsequent relaxation. We intend that the models derived here 
feed a more detailed investigation, with emphasis on the gaseous evolution of the companion 
and a more accurate model for the potential field of Fornax. 

The stellar tails will likely be accompanied by an extensive dispersal of gas; the subse- 
quent infall back to the core occurs over a timescale (Hernquist & Spcrgel 1992, Hibbard & 
Mihos 1995) that far exceeds the time since our postulated merger. In spite of the significant 
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disruption to the disk, the original stellar core of the companion survives in our simulations. 
While a rapid segregation of stars and gas may occur in an isolated galaxy (Weil & Hernquist 
1993), the tidal action on the remnant distribution may make gaseous infall less efficient, 
and thus the shell regions would have exhibited the merger- induced starburst rather than 
the galactic nuclei. Photometric analysis of stars within the inner and outer overdcnsities 
(Coleman et al. 2004, 2005) reveal a small age range centered on 2 Gyr; this truncation of 
star formation may be associated with the quick exhaustion of this cold gas under the action 
of tidal shocks, or instead a rapidly decreasing gas density (Kojima & Noguchi 1997). Where 
infall does occur, the distinct velocity bimodality common to both the simulated Fornax B 
remnant and the metal-poor subpopulation of the Fornax dSph may support the concept of 
a young subpopulation recently introduced to the otherwise aged Fornax dSph. 

Our merger scenario may also add to the discussion on Fornax's globular clusters. 
Firstly, Buonanno et al. (1999) state a dichotomy in age between GCs with H4 being 
uniquely younger than clusters HI, H2, H3 and H5 by ~3 Gyr. This contrasts with the 
opinion of Strader et al. (2003) where H5 is several years the youngest. In either case, the 
age discrepancy may be resolved by the accretion of the younger GC initially in the compan- 
ion galaxy. Based on observations of Forbes and Bridges (2011), several of the Galactic GCs 
have been traced to the Fornax-Leo-Sculptor system, and therefore may have recently been 
part of this companion; the tidal streams we have predicted may contain further GCs as ev- 
idence of this. Secondly, the dynamical heating attributable to a merger event (Kazantzidis 
et al. 2008) could excite the orbits of the GCs. Therefore it is possible that the GCs have 
not since sunken to the core of the host and coalesced into its nucleus due to the past merger 
events of the Fornax dSph. Oh et al. (2000) calculated the timescale for mass accretion to 
the galaxy core (without tidal influences), as ~1 Gyr, therefore implying that this heating 
mechanism may be required at several instances in the past. While this may appear pre- 
sumptuous, it naturally supports the common observation that star formation in Fornax has 
occurred in a sequence of discontinuous starbursts. 

We gratefully acknowledge the suggestions of the anonymous referee that improved this 
manuscript. 
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Fig. 1. — The radial distance between Fornax A and MW (red), Fornax B and MW (blue) for 
P07 (top panel) and Mil (middle panel); the lower panel shows the radial distance between 
Fornax A and B for P07 (Red) and Mil (blue). 



Table 1. Summary of model parameters 



Model 


Proper Motion 


Tm 


Vx 


Vy 


Vz 


rAB 




^s,dw 






Gyr 








kpc 


kpc 


kpc 


1 


P07 


-3.5 


-137 


-80 


103 


3.75 


0.645 


0.6 


2 


P07 


-2.1 


-137 


-80 


103 


3.75 


0.645 


0.6 


3 


P07 


-2.1 


-137 


-80 


103 


3.75 


0.43 


0.4 


4 


Mil 


-3.5 


-226 


-185 


149 


3.75 


0.645 


0.6 


5 


MllL 


-3.5 


-125 


-96 


99 


3.75 


0.645 


0.6 




MllH 




-326 


-274 


199 
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Fig. 2. — The stellar distribution of Fornax B as projected in galactic spherical coordinates 
and centered on the location of Fornax A. The time T (in units of Gyr) is shown at the 
top of each panel. The colour scale in each panel relates to the cumulative number of stars 
within each pixel of the 600x600 resolution image, for Model 1 (top-left, top-center). Model 
2 (lower-left). Model 3 (lower- middle), and Model 4 (lower-right). For comparison, the 
Gaussian-convolved RGB stellar distribution of Fornax (Figure 23 of Coleman et al. 2005) 
is given in the top-right. 
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Fig. 3. — (Left) The stellar (mass) surface density of Fornax B for Models 1, 2 and 3, 
compared with stellar (number) density data (reproduced from Coleman et al. 2005). The 
distribution of stars from each N-body simulation is divided between four on-sky quadrants 
centered on Fornax A, at time T (in units of Gyr). (Right) The star count with Heliocentric 
Velocity for stars within distance bins up to the Fornax tidal radius, for (top) observational 
data reproduced from Battaglia et al. (2006), and (bottom) Model 1 from this work. 
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Fig. 4. — The stellar distribution of Fornax B at present time (T = 0) for Models 1, 3, 4 
and 5. The colour scale is the same as used in Figure 2. The orbit of Fornax A is portrayed 
with the red dashed line. 



